Radiation-reaction force on a particle plunging into a black hole 
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We calculate the self force acting on a scalar particle which is falling radially into a Schwarzschild 
black hole. We treat the particle's self-field as a linear perturbation over the fixed Schwarzschild 
background. The force is calculated by numerically solving the appropriate wave equation for each 
mode of the field in the time domain, calculating its contribution to the self force, and summing 
over all modes using Ori's mode-sum regularization prescription. The radial component of the force 
is attractive at large distances, and becomes repulsive as the particle approaches the event horizon. 



PACS number(s): 04.25-g, 04.70.-s, 04.70.Bw 

The problem of finding the equations of motion for 
a particle in curved spacetime is a long-standing open 
problem in General Relativity. Recently, this problem 
has also become timely and extremely important. The 
planned Laser Interferometer Space Antenna (LISA) is 
expected to detect (among other sources) the gravita- 
tional waves emitted from a compact object orbiting a su- 
permassive black hole (BH). Accurate templates, which 
include also the radiation-reaction (RR) effects on the 
compact-object's orbit, are essential for the detection of 
the signal. 

The traditional approach for calculation of the orbital 
evolution under RR requires the calculation of the fluxes 
at infinity and through the BH's event horizon (EH), of 
quantities which are constants of motion in the absence 
of RR. Then one uses balance arguments to relate these 
fluxes to the local quantities of the object jy]. However, 
such techniques typically fail, because the evolution of 
the Carter constant, which is a non- additive constant of 
motion, cannot be found by balance arguments Q. 

Several prescriptions to include the RR effects in the 
orbital evolution have been suggested. Quinn and Wald 
H and Mino, Sasaki, and Tanaka || recently proposed 
general approaches for the calculation of self forces. How- 
ever, it is not clear how to practically apply these ap- 
proaches for actual computations, the greatest problem 
being the calculation of the non-local "tail" contribution 
to the self force, which arises from the failure of the Huy- 
gens principle in curved spacetime. More recently, Ori 
proposed a practical approach for the calculation of the 
self force which is based on decomposition of the 

self force into modes, and on a mode-sum regularization 
prescription (MSRP). MSRP has been developed in full 
in Refs. [0J|] for a scalar particle in static spherically- 
symmetric spacetimes, and has been applied for several 
non-trivial cases including a static scalar charge outside 
a Schwarzschild black hole (SBH) || and a scalar charge 
in circular orbit around a SBH . In addition, there is 
strong evidence that MSRP is applicable also for electric- 



field RR HQ, and some evidence that it is applicable 
also for gravitational- field RR jl^] . 

MSRP has been directly applied until now only for 
stationary problems, were the field was decomposed into 
Fourier- harmonic modes, and the analysis was done in 
the frequency domain. This was easy to be done in Jl(j] 
for the case of circular orbits around a SBH, because 
the RR in that case admits a discrete spectrum. How- 
ever, in general one faces a time-dependent, evolutionary 
problem, and one expects the spectrum to be continuous 
rather than discrete. In this paper, we apply MSRP for 
the first time to a time-dependent, dynamical problem. 

We consider a pointlike massless particle of scalar 
charge q, moving along a radial (timelike) geodesic out- 
side a SBH of mass M 3> \q\, where the metric is 
ds 2 = -F(r)dt 2 +F- 1 (r)dr 2 + r 2 dft 2 , dtt 2 being the 
metric on the unit 2-sphere, and F(r) = 1 — 2M/r. Let 
the particle's worldline be represented by = x^(t), 
with t being the proper time along the geodesic. For 
inward radial geodesic motion, to be considered here, we 
have (in Schwarzschild coordinates) p — ip p = 0, 



I 1/2 



and i p = E/F(r p ), (1) 



where a dot denotes d/dr, and E is the energy parameter 
(which is a constant of motion in the absence of the self 
force). The scalar field <& coupled to the particle satisfies 
the inhomogeneous wave equation 



6 4 [x^-x^r)](-g)-^ 2 dr, (2) 



g being the metric determinant, and □ denoting the co- 
variant wave operator. We next decompose $ into modes 



£=0 



if 



(3) 



where Yi m {0,Lp) are the standard scalar spherical 
harmonics, an asterisk denotes complex conjugation, 
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and j/r = ip e (r 7 1; r p , t p ). By expanding the delta 
function in Eq. (0) as <5(6* — 6 p )§(tp — tp p ) — 
Y.i,m s ' mOY tm{0,(p)Y l * m {O l>1 if v ) and using the orthogo- 
nality of the Yi m , we find that satisfies 



y uv + V e (r)4> e = S(r;r p ). 



(4) 



Here, w = t + r* and u = t — r* (with r* = r + 2M ln[(r - 
2M)/2M]) are the ingoing and outgoing Eddington null 
coordinates, correspondingly, the effective potential is 
given by V e (r) = (F/4) [£(£ + l)r- 2 + 2Mr~ 3 ] , and 



S=f[ S{r-r p )S(t-t p )dr = 



27e 



*(r-r p ). (5) 



[In the last equality we use dr = dt/t = (F/E)dt, fol- 
lowed by integration over t.] Finally, we express the 
modes <j> 1 in terms of the (m-independent) functions ip e 
by summing over the azimuthal numbers m in Eq. (|^). 
For radial motion we thus find that 



1W 
2 r 



(6) 



The total regularized self force (including both the local 
and the tail parts) exerted on the scalar particle, / RR , 
can be calculated by [0j§O] 



r 

■I a 



RB 



E 



(7) 



(evaluated on the particle's worldline), where /~ = qcj) a 
is the (covariant) self-force contribution associated with 
the ^-mode of the particle's self-field, and and B a 
are regularization parameters, whose values are given by 
Eqs. (101) and (134) of Ref. [§. For the radial geodesies 
considered here, these parameters take the form 



B r ,=-^{2F-E>) 



Af = ±^r, (8) 



(with all quantities evaluated at x^ — xp 1 ). One should 
use fi + ,A+ (fi~,A~) when one calculates the field's 



gradient from the r 



(r 



limit (in general 



fa + + /a" W- (° f course ' the Physical quantity /£ R 
can be derived from either of these two values, or from 
any of their linear combinations.) In practice, we take 
below the r — > r~ limit. 

Thus, in practice, to derive the self force along any 
given radial geodesic (parametrized by E), one should 
first solve Eq. (^) for the various modes (with appro- 
priately chosen initial data — see below), then construct 
the quantities f a , and finally sum over the regularized 
self-force modes using Eq. (@). This sum over modes is 



expected to convergence at least as l/£, as the 0(l/£) 
term in the 1 /£ expansion of f a vanishes . 
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FIG. 1. Numerical grid cell containing a section of the 
worldline. (See the text for further explanation.) 

To solve for tp e , we integrate Eq. (Q) numerically (in the 
time domain) on a double-null grid. This grid is spanned 
by v and u, covering the entire exterior of the SBH (with 
the EH approached at u — > oo). A characteristic initial- 
value problem for ip £ is set up by specifying initial data 
on two null hypersurfaces v = vq and u = uq, taken 
to intersect at some point along the particle's worldline. 
As initial data we take the exact solution corresponding 
to a static particle held fixed at = (vq — uq)/2 ||. 
This solution is not the actual initial field of the geodesic 
particle considered here (which is unknown, in general). 
However, it does approximates the initial field if Tq is 
chosen to be a turning point of the geodesic (when such 
exists) (these initial data are inexact because the accel- 
eration of the geodesic particle at t$ is not the one of 
a static particle, although its position and velocity are), 
or — for a marginally bounded particle (with E = 1) — if 
r-Q is taken large enough. The difference between the ac- 
tual initial field and the static initial data results in the 
occurrence of spurious waves superposed on the actual 
field; however, one may expect such waves to quickly die 
off, unveiling the intrinsic behavior of the field. Numeri- 
cal experiments showed that this is indeed the case: The 
spurious waves were found to decay fast in all cases ex- 
amined (see Fig. ^ e.g., for 7q = 40M). For a marginally 
bounded particle it has been confirmed that the field left 
after the spurious waves decay becomes independent of 
Tq — indicating that one indeed extracts the actual phys- 
ical behavior. (In addition, we found that the larger rg, 
the smaller the amplitude of the spurious waves, and the 
quicker they decay.) 

To construct the difference scheme for the numerical 
integration we use a method similar to that applied by 
Lousto and Price in Ref. fli|] . We integrate the field 
equation (|j) over the unit cell shown in Fig. [l], which 
is centered at v, u and whose sides are of length 2h. 
Let ip\ = il) l (v — h,u — h), ip2 = ip (v — h,u + h), 
ip3 = t//(v + h,u — h), and ^4 = iji e (v + h, u + h), and 
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suppose that ipi, ip2, and 03 are already known, and we 
wish to calculate 04- Integration over the ijj uv term in 
Eq. (0) yields (exactly) ipi — 02 — "03 + 04- Integration 
over the potential term yields (0i + 04) [1 + h 2 V l {r)\ — 
0>2 + ip3)[l - h 2 V l {r)] + 0(h 3 ). (Note that because ip e 
is continuous across the worldline, the integration of the 
potential term here is much simpler than in 113] , where 
the metric perturbations were studied using the Moncrief 
gauge, in which the wave function suffers a discontinuity 
across the worldline.) Finally, integrating over the source 
S (which is most easily done by transforming to the r, t 
coordinates, recalling that dvdu — 2F~ 1 dr dt), we ob- 
tain Z = J S dv du = if the worldline does not cross 
the cell, or, if it does, 

Z = E~ l [k(t ont ) - fc(tin)] (*out - tin) + 0{h 3 ). (10) 

Here, k{t) = F[r p (t)]/r p (t), and t in (t out ) is the t value 
where the worldline enters (leaves) the cell. We can now 
extract the desired quantity 04. To 0(h 2 ) we find 

04 = -01 + [1 - 2h 2 V\r)] (0 2 + 3 ) + Z. (11) 

Our code, which is second-order convergent, evolves the 
scalar field in a straightforward marching. At each grid 
cell, the code does the following: (i) it decides whether 
or not the given worldline crosses the cell; (ii) if the 
worldline crosses the cell, it determines the point where 
it leaves it and calculates t out (given ti n ) to 0(h 2 ); (iii) 
it uses Eq. ( |ll| ) to calculate the field f at the cell's up- 
per point; and (iv) at grid cells containing a section of 
the worldline, it constructs the quantities /~, and fi by 
appropriately extrapolating the field gradients along the 
worldline (based on the already-derived values of the field 
at a few neighboring grid points). 
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FIG. 2. Fhe individual modes of the regularized RR force 
for a particle released from rest from Vq = 40M as a function 
of r/M. Shown are the first 18 modes {£ = 0..17). The 
monopole (£ = 0) modes are displayed by dashed lines, and 
the dipole (I — 1) modes by thick lines. Except for the £ — 0, 1 
modes, the modes' amplitudes decrease monotonically with t. 
Top panel (A): f r ,. Bottom panel (B): ff. 



We next present our results for a particle released from 
rest at Tq = 40M (similar results are obtained also for 
other values of Tq and for the marginally-bound case). 
Figure I displays the behavior of the f r , (§A) and f% (§B) 
components of the RR force. The £ = components are 
everywhere negative, whereas all the other modes (£ > 1) 
are everywhere positive. Figure ^| also shows the decay of 
the spurious waves. Clearly, for values of r smaller than 
25M they are already too small to be noticed. 
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FIG. 3. Top panel (A): The relative difference (/'»///) - 1 
for t = (upper) and I = 1 (lower) as functions of 
(r - 2M)/M. Bottom panel (B): The behavior of the individ- 
ual modes as a function of the mode number £ for different 
values of r. Shown are /*« for r — VIM (*), r — 6M (+), and 
r = 2.01M (x), and ff for r = VIM (O), r = 6M (o), and 
r = 2.01M (□). 

Three properties of the behavior of the individual 
modes are particularly interesting: First, the dipole 
{£ = 1) modes behave differently than the other modes, 
and the closer to the BH, the less important they are. 
Second, the relative importance of the higher modes in- 
creases approaching the BH. This will require care in 
the evaluation of the remainder of the series when we 
sum over all modes (see below). Third, as we approach 
the BH — > f^. This latter property is obvious from 
the following consideration: The covariant components 

and ffj [where U is the outgoing Kruskal coordinate, 
satisfying U oc exp(— u/4M) near the EH] assume finite 
values at the EH itself, as v and U are regular coordi- 
nates at the EH. Consequently, vanishes exponentially 
with u approaching the EH, yielding ff + 
and = f* — /„—>■ fy as we approach the EH. Thus, 
— > //. This is, in fact, a result of spatial gradients be- 
coming comparable to temporal gradients near the EH. 
This behavior is shown in Figure [|(A) for two modes 
(£ = 0, 1), but similar behavior is found also for all the 
other modes. Figure ||(B) displays the behavior of the 
modes as a function of the mode number £, for various 
values of r. The individual modes behave like l~ 2 for 
large values of £. Note, that the closer the particle to the 
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BH, the later the asymptotic l~ 2 behavior starts. Most 
importantly, the detailed behavior of the modes confirms 
the expressions for the analytically-derived regularization 
parameters |7j||. 
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FIG. 4. The full RR force as a function of r. Top panel 
(A): / RR for a marginally-bound worldline (solid line) and 
the curve —0.1 x (r/M)~ 3 (dashed line). Bottom panel (B): 
Free fall from rest starting from = 40M. Dashed line: / t RR . 
Solid line: / R , R 

Next, we sum over all modes to find the total RR 
force. As noted above, the relative importance of the 
higher modes increases as we approach the horizon. This 
causes two problems: (i) it is crucial to include an accu- 
rate approximation of the remainder of the series due to 
our computation of only a finite number of modes, and 
(ii) the noise contribution from the ^-mode to the overall 
force increases with I. The i~ 2 behavior of the modes 
indicates that we can sum over the modes and calculate 
the remainder as was done in Ref. |l(J. Specifically, the 
full RR force is 



/™ = E^ (rcs) +^ 



e+i 



(12) 



where the remainder can be approximated by IZ^ 1 ~ 
^ 2 /" ai V (1) 0?+ !)• Here, f£ (rcg) is the regularized t- 
mode of the force, and il>^(x) w x~ x + x~ 2 /2 + x~ 3 /6 + 
0(x~ 5 ) is the trigamma function. As we sum the series 
only up to £ = 17, this approximation for TZ^ 1 guar- 
antees accuracy of 7 x 10~ 4 (we neglect here the contri- 
bution to IZj^ 1 from terms which scale like £~ 3 ). Ob- 
viously, approaching the EH / r » — > f t (as each of the 
individual modes does). Figure^ shows the full RR force 
as a function of r for two cases: Fig. |(A) shows / R . R 
for a marginally-bound trajectory (E = 1). At large dis- 
tances this force behaves like / RR « —(G/c 2 )(3q 2 M/r 3 . 
The exponent of r is found here to a 1% accuracy, and 
we find the parameter (3 — (1.00 ± 0.15) x 1Q -1 . Fig- 
ure |](B) shows the case of fall from rest, starting from 
Tq = 40Af. At large values of r both components of the 



force vanish, in accord with the vanishing of the force for 
a static scalar charge. The covariant t component, / RR , 
is everywhere positive and increases monotonically ap- 
proaching the BH. This is a consequence of the particle 
losing energy by radiating, part of which escapes to infin- 
ity, and the rest being captured by the BH. The covari- 
ant r* component, / RR , is attractive at large distances. 
However, near the peak of the effective potential barrier 
(near r* « 0) its behavior changes, and near the EH 
it approaches the value of / RR , as expected. Note that 
both components arrive at the EH at a bounded value. 
Because / t RR is expected to be positive (the particle only 
loses energy by radiating), we infer that / RR would also 
be positive approaching the EH, under very general con- 
ditions. In particular, if this behavior persists also for 
charged BHs, and for an electrically-charged particle, 
then the properly-defined covariant spatial-component of 
the RR force at the EH would be repulsive. If this is 
indeed the case, then the RR force acts to reduce the 
parameter space for which a nearly-extreme spherical 
charged BH can be overcharged [jl5). The question of 
whether cosmic censorship for that case is saved by RR 
effects, however, awaits further considerations. 
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